Representation of conformal maps by rational functions

The traditional view in numerical conformal mapping is that once the boundary correspondence function has been found, the map and its inverse can be evaluated by contour integrals. We propose that it is much simpler, and 10-1000 times faster, to represent the maps by rational functions computed by the AAA algorithm. To justify this claim, first we prove a theorem establishing root-exponential convergence of rational approximations near corners in a conformal map, generalizing a result of D. J. Newman in 1964. This leads to the new algorithm for approximating conformal maps of polygons. Then we turn to smooth domains and prove a sequence of four theorems establishing that in any conformal map of the unit circle onto a region with a long and slender part, there must be a singularity or loss of univalence exponentially close to the boundary, and polynomial approximations cannot be accurate unless of exponentially high degree. This motivates the application of the new algorithm to smooth domains, where it is again found to be highly effective

New Laplace and Helmholtz solvers

New numerical algorithms based on rational functions are introduced that can solve certain Laplace and Helmholtz problems on two-dimensional domains with corners faster and more accurately than the standard methods of finite elements and integral equations. The new algorithms point to a reconsideration of the assumptions underlying existing numerical analysis for partial differential equations

Efficient algorithms for computing rank-revealing factorizations on a GPU

Standard rank-revealing factorizations such as the singular value decomposition and column pivoted QR factorization are challenging to implement efficiently on a GPU. A major difficulty in this regard is the inability of standard algorithms to cast most operations in terms of the Level-3 BLAS. This paper presents two alternative algorithms for computing a rank-revealing factorization of the form $A = U T V^*$, where $U$ and $V$ are orthogonal and $T$ is triangular. Both algorithms use randomized projection techniques to cast most of the flops in terms of matrix-matrix multiplication, which is exceptionally efficient on the GPU. Numerical experiments illustrate that these algorithms achieve an order of magnitude acceleration over finely tuned GPU implementations of the SVD while providing low-rank approximation errors close to that of the SVD

Fullwave design of cm-scale cylindrical metasurfaces via fast direct solvers

Large-scale metasurfaces promise nanophotonic performance improvements to macroscopic optics functionality, for applications from imaging to analog computing. Yet the size scale mismatch of centimeter-scale chips versus micron-scale wavelengths prohibits use of conventional full-wave simulation techniques, and has necessitated dramatic approximations. Here, we show that tailoring "fast direct" integral-equation simulation techniques to the form factor of metasurfaces offers the possibility for accurate and efficient full-wave, large-scale metasurface simulations. For cylindrical (two-dimensional) metasurfaces, we demonstrate accurate simulations whose solution time scales \emph{linearly} with the metasurface diameter. Moreover, the solver stores compressed information about the simulation domain that is reusable over many design iterations. We demonstrate the capabilities of our solver through two designs: first, a high-efficiency, high-numerical-aperture metalens that is 20,000 wavelengths in diameter. Second, a high-efficiency, large-beam-width grating coupler. The latter corresponds to millimeter-scale beam design at standard telecommunications wavelengths, while the former, at a visible wavelength of 500 nm, corresponds to a design diameter of 1 cm, created through full simulations of Maxwell's equations.Comment: 11 pages, 6 figure

High-order numerical methods for scattering problems

Scattering problems are ubiquitous in scientific and engineering applications. These problems can often be straightforwardly modeled using relatively simple linear partial differential equations, such as the Helmholtz equation. The numerical solution of these equations, however, can be challenging due to several features of the problem, such as ill-conditioning, unbounded domains, corner singularities, and oscillatory solutions. In this thesis we present novel numerical methods for the efficient handling of some of these challenges. We first consider the case of homogeneous scattering on domains with corners and introduce a new solver in this regime, based on results in the rational approximation theory literature. The solver uses an exponentially clustered sum of dipoles to resolve the singularities in the solution and a simple expansion to resolve the smooth part. We show that the convergence is root-exponential with respect to the number of degrees of freedom and illustrate the behavior of the solver in several numerical experiments. We next turn our attention to the case of inhomogeneous scattering and present a new solver for the Lippmann–Schwinger equation, an integral equation reformulation of the problem. The new solver is based on the hierarchically block separable matrix format and exploits the geometry of the discretization and problem for accelerated inversion. The solver is shown to be effective both as a direct solver and preconditioner in a number of numerical experiments. For our last contribution we consider broadband applications in homogeneous scattering. We present a new technique for accelerating existing direct solvers in this regime. The technique works by computing basis matrices for all of the low-rank blocks in a rank-structured matrix format and leverages the observation that the resulting ranks are often comparable to those at just the highest frequency. Again, the viability of the method is illustrated in several numerical experiments